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Abstract 

Purpose Despite the fundamental role of ecosystem goods 
and services in sustaining human activities, there is no 
harmonized and internationally agreed method for including 
them in life cycle assessment (LCA). The main goal of this 
study was to develop a globally applicable and spatially 
resolved method for assessing land use impacts on the 
erosion regulation ecosystem service. 

Methods Soil erosion depends much on location. Thus, unlike 
conventional LCA, the endpoint method was regionalized at 
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the grid cell level (5 arcmin, approximately 10 x 10 km ) to 
reflect the spatial conditions of the site. Spatially explicit 
characterization factors were not further aggregated at broader 
spatial scales. 

Results and discussion Life cycle inventory data of topsoil 
and topsoil organic carbon (SOC) losses were interpreted at 
the endpoint level in terms of the ultimate damage to soil 
resources and ecosystem quality. Human health damages 
were excluded from the assessment. The method was tested 
on a case study of five 3-year agricultural rotations, two of 
them with energy crops, grown in several locations in Spain. 
A large variation in soil and SOC losses was recorded in the 
inventory step, depending on climatic and edaphic condi¬ 
tions. The importance of using a spatially explicit model and 
characterization factors is shown in the case study. 
Conclusions The regionalized assessment takes into ac¬ 
count the differences in soil erosion-related environmental 
impacts caused by the great variability of soils. Taking this 
regionalized framework as the starting point, further re¬ 
search should focus on testing the applicability of the meth¬ 
od through the complete life cycle of a product and on 
determining an appropriate spatial scale at which to aggre¬ 
gate characterization factors in order to deal with data gaps 
on the location of processes, especially in the background 
system. Additional research should also focus on improving 
the reliability of the method by quantifying and, insofar as it 
is possible, reducing uncertainty. 

Keywords Ecosystem services • Land use impacts • 
Regionalized life cycle impact assessment • Soil organic 
carbon • Soil loss 

1 Introduction 

Life cycle assessment (LCA) aims to provide a general picture 
of the environmental impacts of resource consumption and 
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emissions during the entire cycle of a product system. LCA 
normally uses spatially and temporally independent linear 
impact assessment (LCIA) models to calculate the environ¬ 
mental damage caused by the resource consumption and 
emissions quantified in the life cycle inventory (LCI) of the 
product system under analysis. Site-specific conditions are 
therefore not usually taken into account when the environ¬ 
mental impacts of the product are evaluated. The impacts 
of the different LCA categories have consequences for the 
environment and human welfare on different spatial scales. 
The spatial scale at which these consequences occur has 
nothing to do with the importance of the categories, but 
with a need for spatial differentiation for some impact 
categories. Since economic processes are spread world¬ 
wide, local impacts have a global extension as well 
(UNEP 2003). 

Despite the fundamental role of ecosystem goods and 
services in sustaining human activities, such as climate regu¬ 
lation, nutrient cycling, and erosion control (Millennium 
Ecosystem Assessment 2005), there is no harmonized and 
internationally agreed method for including them in LCA. 
According to a key framework for land use impact assessment 
in LCA (Mila i Canals et al. 2007a), ecosystem goods and 
services that should be integrated within LCA are impacts on 
biodiversity and, at least, impacts on the following five major 
ecosystem services: biotic production potential, carbon se¬ 
questration potential, freshwater regulation potential, water 
purification potential, and erosion regulation potential. 
Operational characterization factors and methods have been 
developed for impacts on biodiversity (De Schryver et al. 
2010; Koellner and Scholz 2008) and on ecosystem services 
(Brandao and Mila i Canals 2012; Muller-Wenk and Brandao 
2010; Saad et al. 2011). The methodology for assessing the 
degradation of the erosion regulation function due to land 
transformation and land occupation was developed using dif¬ 
ferent spatial scales (Saad et al. 2011). However, erosion 
regulation characterization factors (CLs) were developed 
within a Canadian context and therefore are not globally 
applicable (Saad et al. 2011). 

Because an activity’s land use impacts depend on local 
conditions, a conventional site-independent LCA methodol¬ 
ogy might not be accurate. Methods focusing on these 
impacts should therefore include geospatial information in 
both the LCI and LCIA phases. Different levels of region¬ 
alization (e.g., countries, ecoregions, and biomes) and eco¬ 
logical unit classifications (e.g., life zones by Holdridge 
1947 and ecoregions by Olson et al. 2001) in the LCIA of 
a range of land use impacts are presently used without a 
clear recommendation on a standardized approach to ad¬ 
dress spatial differentiation. Although the use of ecological 
or geographical units instead of administrative borders pro¬ 
vides better estimates of the site dependency of land use 
impacts, especially in countries with a high degree of 


variability, it is generally easier to find information to create 
regionalized impact assessment CFs at the country scale. 

Our objective was to go one step further toward the inte¬ 
gration of ecosystem services in LCA by developing a glob¬ 
ally applicable and spatially resolved method to include land 
occupation impacts on the erosion regulation. Indicators of the 
impact category were defined on the endpoint level and were 
modeled up to the entities described by the areas of protection 
(AoP), i.e., natural resources, natural environment, and human 
health. The case study conducted to demonstrate the applica¬ 
bility of the method focused on the impacts of agricultural 
rotations with energy crops in Spain as compared to the 
cultivation of rotations with only traditional crops. 

2 Materials and methods 

Figure 1 illustrates the general land use impact mechanism 
and shows the three impact pathways studied as follows: 

Land occupation leads to soil erosion, and this leads to loss 
of topsoil reserves, which leads to soil resource depletion 
(AoP natural resources, impact pathway 1 in Fig. 1). 

Land occupation leads to soil erosion and altered soil 
function, which affects net primary production and 
leads to damage to ecosystem quality (AoP natural 
environment, impact pathway 2 in Fig. 1). 

Land occupation leads to soil erosion and altered soil 
function, which affects net primary production and 
leads to damage to human health (AoP human health, 
impact pathway 3 in Fig. 1). 

2.1 Soil resource depletion impact pathway 

This AoP represents damage from the impact of removal of 
soil resources from the environment on soil resource deple¬ 
tion (impact pathway 1 in Fig. 1). Annually, humans cause 
the loss of 50-75 billion metric tons of soil (Harvey and 
Pimentel 1996). Agricultural land accounts for 75 % of the 
soil erosion worldwide, though it also occurs in other human- 
modified ecosystems, such as during the constmction of roads 
and buildings. More than 75 % of the arable soils of the world 
suffer from moderate to very high soil losses (Reich et al. 
2001), generally ranging from 10 to 100 tha _1 y ear_1 
(Pimentel et al. 1987), which is at least tenfold higher than 
the average rate of soil formation of 0.5-1 tha'year -1 (Mann 
et al. 2002). Current soil losses due to land use (land occupa¬ 
tion) reduce soil availability as a future resource. 

2.1.1 Life cycle inventory data requirements 

The type of land use has a determining role in the quantity 
of soil loss as specific direct physical interventions are often 
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a Land occupation does not “cause changes” but contributes to prolong altered conditions 

b Benefits from downstream deposition of eroded material disregarded in the paper, as sediments cross 
the system boundary 

Fig 1 Main impact pathways related to land use. The pathways discussed in the paper are shown with cross-hatching and thicker arrows (adapted 
from the ILCD handbook; JRC 2010) 


related to land management (Frischknecht et al. 2007). In 
LCA, the ecoinvent database v3.0 (Weidema et al. 2011), 
based on Global Land Cover 2000 (Bartholome and 
Belward 2005) and ecoinvent v2.0 (Frischknecht et al. 
2007), identifies ten main types of land use and land cover 
classes in the first level (e.g., agriculture, forest) and pro¬ 
vides more detailed information on land intensity and man¬ 
agement in the following levels (e.g., arable non-irrigated 
and arable non-irrigated extensive). This tiered structure 
allows for different levels of detail in the LCI, depending 
on the quality requirements of the LCA study and the 
information available on the product under study. We chose 
this flexible classification system to record the land occupa¬ 
tion type in the LCI. Yet, ecoinvent v3.0 does not distinguish 
between cultivated crops and specific crop management prac¬ 
tices (e.g., type and timing of tillage operations and use of 
residue mulches), which are key factors controlling soil ero¬ 
sion on agricultural lands. 

Soil losses due to land occupation must also be included 
in the LCI as soil loss mass (in grams). There are many 
estimation models that can be used to estimate soil loss, 
such as PESERA (Gobin and Govers 2003), INRA (Le 


Bissonais et al. 2002), or universal soil loss equation erosion 
model (USLE; Wischmeier and Smith 1978). From all of 
them, we recommended the USLE to register soil losses in 
the LCI. USLE takes into account the effect of a particular 
land use type on soil erosion by water. There is a consensus 
that the USLE equation and its update (revised USLE, 
RUSLE; Renard et al. 1997) are valid methods to estimate 
soil losses by water at the inventory stage (Beck et al. 2010; 
Muys and Garcia Quijano 2002). The use of (R)USLE to 
predict soil erosion losses has increased worldwide due to 
the growing availability of and accessibility to climatic, 
edaphologic, and land use and land cover data at the local 
and regional levels. This spatial information can be increas¬ 
ingly found in a geo-referenced format, thus allowing for 
data processing and visualization in a geographic informa¬ 
tion system (GIS) software. Furthermore, there are many 
biophysical models, such as EPIC (Williams et al. 1984) or 
APEX (Williams and Izaurralde 2005), which already have 
algorithms to simulate soil erosion with (R)USLE for many 
crops and at different spatial resolutions. Results from sce¬ 
nario simulations from these models can be incorporated in 
the inventory stage to register soil erosion. 
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Soil formation was not contrasted to soil erosion in this 
study because factors other than human land use, namely, 
climate, and the soil parent material, are generally recog¬ 
nized as the most important aspects governing soil forma¬ 
tion (Jenny 1994). Data on the area as well as the duration 
required for the production of a certain amount of products 
and services have to be collected in the LCI as land occu¬ 
pation impacts are recorded in area per time (e.g., in square 
meters per year). The geo-referenced location of land use 
(longitude/latitude) should also be included if available. 
Failing that, a broader resolution (e.g., region, country) 
can be used, although this reduces the quality of the LCI 
data by increasing uncertainty due to the high variability of 
spatial conditions. 

In summary, the following elementary flows need to be 
accounted in LCI: type of land use, soil erosion by water, time 
of occupation, area of occupation, and location of the activity. 

2.1.2 Impact assessment model 

Similar to other endpoint methods (Ecoindicator 99, EI99; 
Goedkoop and Spriensma 2001), damage to resources is 
expressed as surplus energy needed to make the resource 
available at some point in the future. This is a suitable unit 
to evaluate soil depletion, which indicates the anticipated 
energy removal from nature to provide a unit of soil eroded 
during land occupation. Instead of using energy units (meg¬ 
ajoule equivalents), such as in EI99, we used emergy units 
(megajoule solar equivalents). Emergy is the environmental 
work needed for the formation of a natural resource (Odum 
1996). In emergy, solar energy is the unit of reference as it is 
the primary energy source on the earth. Unlike the energy 
metric, emergy accounts for quality differences of the ener¬ 
gy used to generate a product or service by converting raw 
units (e.g., kilograms soil and cubic meter water) to a 
common basis, i.e., units of solar energy. For example, 
sunlight, soil, water, and coal are assessed together by 
expressing them all with the solar energy unit. The advan¬ 
tage of emergy compared with other energy and exergy units 
is that it not only accounts for energy carriers (e.g., gas) and 
non-energetic materials (e.g., minerals) but also assigns an 
exergy (useful energy) value to land use (Rugani et al. 
2011). Therefore, emergy evaluates the (solar) energy the 
natural system is deprived of to yield the new stock of soil 
lost during land use. We assigned an average transformity 
(i.e., the emergy required to make available one unit of 
energy, mass, or volume of the resource) for all soil types 
and locations: 23.9 MJ se g _1 soil (Odum 1996). This means 
that around 24 MJ of solar energy is required to generate a 
gram of soil lost by erosion. In LCA, transformity has 
been called solar energy factor (SEF; Rugani et al. 2011). 
The effect of soil erosion on soil resource depletion (A R) 
is expressed as follows with units of megajoule solar 


equivalent soil loss per unit of area and time of land 
occupation: 

If soil loss = 0, 

If soil loss > 0, 


This endpoint indicator combines the inventory flow (i.e., 
soil loss) with the local available soil reserves (SD 7 , soil 
depth in the specific location i) on a spatial resolution of 
5 arcmin (approximately 10 x 10 km 2 ; FAO/UNESCO 2007) 
and with the solar energy factor of soil (SEF soil ) as the 
characterization factor. Soil loss mass of the LCI is weighted 
with soil depth as the environmental significance of soil loss 
depends on the soil stock size at the specific location. 
Twenty-one soil depth classes were distinguished based on 
the FAO’s soil depth classification map, from very shallow 
soils (0.05 m deep) to very deep soils (2.25 m deep). These 
classes resulted from the combination of five major soil 
depth categories (further details in Electronic supplementary 
material (ESM), Section 1). Local soil reserves on each grid 
cell of the soil depth map were normalized with a reference 
soil depth (SD ref ). We selected the upper limit of the very 
deep soil depth category as the reference (SD ref =3 m) as an 
indicator of the potential soil quality. This would only occur 
if the whole of the grid cell area had this maximum soil 
depth, which we assumed was not possible once the first soil 
particle was eroded. Note, therefore, that SD 7 can never be 
equal to SD ref , so the characterization factor always takes 
positive values and any amount of soil loss will have an 
environmental impact. Choosing the maximum soil depth as 
the reference was judged to be a more representative cor¬ 
rection of the site dependency of the characterization factor 
rather than the soil reserves of a particular region (e.g., 
Swiss lowlands for ecosystem quality assessment in the 
EI99 methodology; Goedkoop and Spriensma 2001). 
Larger damage factors were assigned to thinner soils by 
calculating the difference between SD ref and SD 7 before 
normalization with the soil reference. This means that, based 
on the same amount of soil loss, thinner soils are more 
vulnerable than thicker soils. The indicator relates impact 
assessment to biogeographical conditions in each grid cell i 
(SD 7 ) without any further aggregation of land use type or 
land use cover. For a land use activity under study, the lower 
the indicator result, the less soil resource depletion and the 
less damage to the environment. 

2.2 Ecosystem quality impact pathway 

This AoP represents damage from the negative impacts of 
altered soil function and net primary productivity on the 
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function and structure of natural ecosystems (impact pathway 
2 in Fig. 1). In the erosion process, soil quality declines, 
essential plant nutrients are lost, and soil depth is reduced. 
As a result, biomass productivity diminishes. Ultimately, this 
can adversely affect overall biodiversity and ecosystem qual¬ 
ity. Numerous positive correlations between plant biomass 
productivity (NPP) and vascular plant species diversity and 
richness have been established (Costanza et al. 2007; 
Flombaum and Sala 2008), though overall biodiversity does 
not always correlate with productivity (e.g., Mediterranean 
hot spots and intensively managed agricultural lands). 

One of the methods used in LCA for measuring environ¬ 
mental impacts on ecosystem functions is based on the soil 
organic matter content (Mila i Canals et al. 2007b). We 
assessed the effects of soil erosion on the terrestrial ecosys¬ 
tem quality by, in a first step (damage factor), linking soil 
loss to soil organic carbon (SOC) loss and, in a second step 
(effect factor), linking SOC loss to biomass production drop. 
Thus, ecosystem biomass production was modeled as a 
function of soil quality, which is indicated by the SOC 
content of the soil lost. We assessed ecosystem biomass 
production as a function of the net primary production of 
potential natural vegetation (NPP 0 , i.e., the anticipated state 
of mature vegetation in the absence of human intervention). 

2.2.1 Life cycle inventory data requirements 

SOC losses must be registered in the LCI. Like soil erosion, 
SOC losses can be directly derived using biophysical mod¬ 
els like EPIC (Williams et al. 1984) or APEX (Williams and 
Izaurralde 2005). Alternatively, SOC losses may be calcu¬ 
lated by, in a first step, determining the topsoil OC content 
and, in a second step, multiplying the already estimated soil 
losses by the percentage of topsoil OC content. The quantity 
of SOC in the topsoil can be determined by direct measure¬ 
ments. In the absence of plot-level soil data, SOC content 
can be determined using already existing spatial data layers 
(e.g., map of organic carbon in topsoils in Europe: Jones et 
al. 2005; Harmonized World Soil Database, HWSD: FAO et 
al. 2009), though inventory data quality will be affected. On 
the other hand, while site-specific measurements are more 
accurate, this is not likely to be achievable in most LCA 
studies (Mila i Canals et al. 2007b). 

The soil unit where the impact activity occurs should be 
included in the LCI using the most recent FAO classification 
(FAO et al. 1990), which identifies 28 different soil units 
(soil types), each with harmonized soil parameters. The soil 
unit can be identified by direct measurements on the occu¬ 
pied land under analysis. Otherwise, it may be approximated 
using literature or spatial databases such as the Harmonized 
World Soil Database (FAO et al. 2009). Again, direct meas¬ 
urements are more accurate, although not always available 
in LCA studies. 


To estimate the SOC in each soil unit, we determined the 
content of the 28 soil units in the FAO classification system 
(Table 1). We averaged the topsoil OC content of over 
16,000 soil mapping units covering the entire land area of 
the world with a resolution of 30 arcsec (approximately 1 x 
1 km 2 ) in the HWSD (FAO et al. 2009). 

The 28 soil units were grouped into seven major catego¬ 
ries (see Table 1). In 24 out of 28 soil units, there was <2 % 
SOC, the threshold selected by the European Commission 
for defining soils in phase of pre-desertification (COM 

2002) . As in the case of the inventory for the soil 
resource depletion indicator, information on the type 
and intensity of land use, area size (in square meters), 
duration (in years), and the location of the occupation 
should be recorded. 

2.2.2 Impact assessment model 

According to the International Reference Life Cycle Data 
System (ILCD) Handbook (JRC 2010), which is a series of 
technical documents that provide detailed guidance on all 
the steps required to conduct a LCA study, species diversity 
is the recommended indicator to be implemented in end¬ 
point LCA methodologies when modeling damage to eco¬ 
system quality. Effects on species diversity are usually 
quantified in terms of the potentially disappeared fraction 
(PDF) of vascular plant species (Koellner 2000), such as in 
the EI99 (Goedkoop and Spriensma 2001), ReCiPe 2008 
(Goedkoop et al. 2009), and IMPACT2002+ (Jolliet et al. 

2003) methods. However, according to ILCD, function- 
related parameters, such as the biomass production of the 
ecosystem, might also be good endpoint indicators. 

A limited number of studies have focused on accounting 
for current or potential NPP losses caused by soil erosion 
due to the complex connection between the two factors. In 
addition to soil properties, NPP depends on many other 
physical environmental aspects, such as leaf area index, 
precipitation, atmospheric C0 2 concentration, and temper¬ 
ature (Melillo et al. 1993). NPP is therefore commonly 
used as an indicator to reflect ecosystem response to 
climate change. Change in climate may decrease NPP 
(lower precipitation or cloudiness) or increase it (photo¬ 
synthesis enhancement). Due to the many relationships 
between factors, estimates of the direct cause-effect chain 
between NPP and soil erosion are scarce and highly 
uncertain. 

Of the studies performed on the relationship between soil 
loss and NPP/NPP 0 loss, most estimate productivity losses 
according to qualitative degrees of erosion (e.g., light and 
extreme) and limit the scope to the local or the regional level 
(Mann et al. 2002; Mokma and Sietz 1992). We developed 
an initial approach to convert soil loss-NPP 0 loss into 
approximate quantitative linear relationships using data at 
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Table 1 Topsoil organic carbon 
(percent weight) of the 28 soil 
units in the Soil Map of the 
World (FAO et al. 1990), 
according to the HWSD (FAO et 
al. 2009), and linear equations 
used in the impact model 


Soil units within the same soil 
category are arranged by in¬ 
creasing SOC content. SOC loss 
in the %NPP 0 depletion (NPPD) 
equations is expressed as grams 
of C loss (in a square meter and 
year) 


Soil units Topsoil %NPP 0 depletion equations 

HWSD organic 

carbon 
(% weight) 


Gypsisols 

0 to <0.5 

If SOC loss = 0, 

NPPD = 

Arenosols 


If SOC loss > 0, 

NPPD = 

Calcisols 



NPPD = 

Solonchaks 




Lixisols 

0.5 to <1.0 

If SOC loss = 0, 

NPPD = 

Luvisols 


If SOC, oss > o, 

NPPD = 

Solonetz 



NPPD = 

Plinthosols 

Planosols 

Fluvisols 

Regosols 

Leptosols 

Acrisols 

1.0 to <1.5 

If SOQoss = o, 

NPPD = 

Vertisols 


If SOQoss > 0, 

NPPD = 

Cambisols 



NPPD = 

Anthrosols 




Kastanozems 




Ferralsols 




Greyzems 

1.5 to <2.0 

If SOQoss = o, 

NPPD = 

Podzoluvisols 


If SOC l 0 s S > 0, 

NPPD = 

Alisols 



NPPD = 

Nitisols 




Phaeozems 




Chernozems 




Gleysols 

2.0 to <2.5 

If SOQoss = 0, 

NPPD = 

Podzols 


If SOC l 0 s S > 0, 

NPPD = 

NPPD = 

Andosols 

4.86 

If SOQoss = 0, 

NPPD = 



If SOQoss > 0, 

NPPD = 

NPPD = 

Histosols 

34.60 

If SOQoss = 0, 

NPPD = 



If SOC l 0 s S > 0, 

NPPD = 

NPPD - 


0 

4.09 x SOQoss + 2.66for SOCi oss < 23.80g 
100 for SOQoss > 23.80g 

0 

1.96 x SOCioss + 2.66 for SOCioss < 49.66 g 
100 for SOQoss > 49.66 g 


0 

1.32 x SOQoss + 2.66 for SOCioss < 73.74 g 
100 for SOQoss > 73.74 g 


0 

0.88 x SOCioss + 2.66 for SOCioss < 110.61 g 
100 for SOQoss > 110.61 g 


0 

0.69 x SOCioss + 2.66 for SOCioss < 141.07 g 
100 for SOQoss > 141.07 g 

0 

0.31 x SOQoss + 2.66 for SOCioss < 314.00 g 
100 for SOQoss > 314.00 g 

0 

0.04 x SOQoss + 2.66 for SOCioss < 2433.50 g 
100 for SOQoss > 2433.50 g 


the global level found in the literature (Dregne and Chou 
1992; FAO/UNEP 1984; Zika and Erb 2009). This was 
done for each soil unit group in Table 1 (see the proce¬ 
dure in the ESM, Section 2). Using this method, soil loss 
was first related to the loss of SOC as a measure of soil 
quality and, finally, to ecosystems’ loss of biomass 
productivity. 

In the impact assessment model, the inventory flow (i.e., 
SOC losses) should be used as an input parameter in the 
equation in Table 1 for the soil unit where the land use 
activity is developed. Note that for SOC losses greater 
than or equal to a predetermined threshold for each soil 
unit group (when NPPD =100), NPP 0 is completely lost 
and the soil is unlikely to be able to recover. Soils with 


low SOC content, such as those in arid and semi-arid 
areas, are less resilient than soils rich in SOC, which are 
usually found in wet regions. These SOC loss thresholds 
are equivalent to a soil depletion of approximately 65 t 
ha -1 y ear_1 > which is extreme soil loss. Most agricultural 
land in the world loses soil at a rate of between 13 and 
40 tha _1 year _1 (Pimentel and Kounang 1998), whereas 
losses of more than 100 tha^year -1 only occur in ex¬ 
treme events (Morgan 1992). 

The effects of soil erosion on ecosystem quality (AEQ) 
are expressed using a growth-based value: NPPD (potential 
net primary production depletion). For an occupation of 
1 m and 1 year, AEQ ranges from 0 to 1 (percentage 
expressed as a decimal). 
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If SOC loss = 0, 
If SOC loss > 0, 


AEQ 

AEQ 


0 

aSOCi oss + b NPPo % 
Axtx-—-x 


100 


NPP 


—v-— 

IX'I 


0,ref 


CF 


NPPD m 2 year 


( 2 ) 


This endpoint model combines the inventory flow (a* 
SOCi 0 ss + b, if SOCi oss >0, i.e., mass of SOC losses trans¬ 
formed into %NPP 0 losses) with NPP 0 values spatially 
resolved for each grid cell (NPPo /) at 5 arcmin (approximately 
10x10 km ; Haberl et al. 2007) to obtain the absolute biomass 
productivity drop at the specific location. Values were then 
normalized with an NPP 0 value corresponding to that of the 
ecosystem with the highest biotic productivity worldwide 
(NPP 0 r ef = 1,496 gCm year ) as a representative reference 
of the potentiality of the ecosystem. Larger impact factors 
were allocated to the most productive soils (higher NPP 0 ,/). 
From an ecosystem services viewpoint, we made the general 
assumption that ecosystems with higher NPP values are con¬ 
sidered more valuable as NPP is a scarce resource on earth 
(Pfister et al. 2011). The most productive lands can also be 
used for a greater diversity of purposes. Unlike Eq. 1, NPP 0 ,/ is 
normalized without first calculating the difference between 
NPP 0 ,ref ^nd NPP 0 ,/. This implies that the CFs of both AoP are 
not at all comparable and do not have to be as they represent 
different environmental concerns. Note that, for complete 
losses of NPP 0 (< 2 x SOC loss +£>= 100 ), damages only depend 
on the NPP 0 ,//NPP 0 ,ref ratio. As for the resource depletion 
indicator, regionalization was at the grid cell level, without 
aggregating values on broader scales. For a land use activity 
under study, the lower the indicator result, the less of a drop in 
biomass production and the less damage to the environment. 


2.3 Case study on energy crops grown in Spain 


The proposed method was applied to agricultural plots of 
food and energy crops in Spain. This was selected as a 
representative case study because soil loss by water erosion 
is one of the main causes of land degradation in the country 
(EEA 2005) and, important for validating the methodology, 
because soil erosion is unevenly distributed, with a gradient 
from north to south and from the Atlantic to the 
Mediterranean coast (MERMA 2012). Because an increase 
in energy crop production is forecast in Spain (EEA 2006), 
it is important to assess the environmental cost of this 
alternative, given the effect of agriculture on the environ¬ 
ment. To illustrate the outlined method, we analyzed the 
environmental impacts of growing five 3-year agricultural 
rotations, two with energy crops, on the erosion regulation 
ecosystem service. The analysis includes 120 agricultural 
plots covering the main Spanish watersheds (ESM Fig. 2SI). 
The functional unit selected for the assessment was a land 


occupation of 1 m in 1 year (square meters per year), con¬ 
sidering the average soil and SOC loss of a complete crop 
rotation system. We used the attributional framework so that 
any potential additional consequence of planting the rotations 
on other parts of the economy was not included, nor did we 
include past impacts on soils from historical land use. Further 
details on the study area, crop production systems, and data 
sources for the LCI are described in ESM, Sections 3 to 5. 


3 Results 

3.1 LCI A characterization factors 

Regional characterization factors for the soil resource deple¬ 
tion and ecosystem quality impact pathways of the erosion 
regulation ecosystem service are shown in Fig. 2. GIS layers 
are also provided in the ESM. For soil resource depletion, 
the ratio is ((SD ref -SD/)/SD ref )x SEF soil . Lower soil depths 
and therefore higher damage factors are found at high northern 
latitudes and over wide areas of Asia. For ecosystem quality, 
Fig. 2 shows the NPPo,//NPPo, re f ratio. Higher biomass pro¬ 
ductivities and therefore damage factors are found at the low 
latitudes of the tropics. For both the soil resource depletion 
and ecosystem quality indicators, lower CF values indicate 
less sensitivity of the ecosystem to potential land use impacts. 

3.2 Case study on energy crops grown in Spain 

Both soil and organic carbon losses and impact factors vary 
as a function of location, thus leading to considerable differ¬ 
ences in the environmental damage from soil erosion in 
different watersheds. 

In the LCI, the crop rotation with the greatest erosion rate 
(Table 2) was when the field lay fallow in the last year 
(winter barley-winter wheat-unseeded fallow). These soil 
losses were around ten times higher than the crop rotation 
with the lowest erosion rate, the energy profitable poplar 
with a short forestry rotation (poplar-poplar-poplar). For 
annual cereal (winter barley-winter wheat-rye), legume 
(winter barley-winter wheat-pea), and energy crop (winter 
barley-winter wheat-oilseed rape) rotation systems, similar 
soil losses were recorded, with rates about 40 % lower than 
for fallow rotation. Outcomes of the LCIA phase showed the 
same pattern as in the LCI phase, as shown in Table 2. 

According to the presented results, the watersheds in 
Spain with the lowest water erosion rates and environmental 
damage, thus making them the most appropriate for rotating 
poplar and oilseed rape energy crops, are the Duero 
(Northern Spain) and the Guadiana (Central Spain), while 
those with the highest water erosion rates and environmental 
impact are the internal watersheds of Catalonia (Northeast 
Spain) and the Jucar watersheds (Eastern Spain). 
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a 

6 Soil resource depletion [MJ se g 1 soil loss] 24 



b 

0 Ecosystem quality [-] 


1 




Fig. 2 Characterization factors for soil resource depletion (a) and ecosystem quality (b) 


4 Discussion 

4.1 Case study on energy crops grown in Spain 

The results of the case study showed that the implementa¬ 
tion of the poplar energy crop rotation system in Spain can 
potentially reduce erosion rates and related environmental 
impacts per area-time unit compared to traditional cereal and 
legume crop rotations in the country. In relation to the 
recommended locations to cultivate the analyzed crop rota¬ 
tions, the results indicate a trend. However, these results 
should be interpreted with caution. They are valid for the 
analyzed plots, but each specific case should be studied sep¬ 
arately as soil erosion depends on plot-level factors. Apart 
from erosion, other aspects such as biodiversity impacts or 


replacement of other crops, which then need to be imported, 
would need to be considered to assess if planting poplar in 
Spain is an environmentally preferable decision. 

Although the area-time unit does not relate to energy crops’ 
function of producing bio-based electricity or fuels, the frame¬ 
work is still meaningful when ranking crop rotations by their 
impact on soil replacement energy and lost productivity. 
Sustainable land management practices reduce impact inten¬ 
sity in a given area. However, intensive agricultural produc¬ 
tion can diminish soil quality and jeopardize the preservation 
of productive agricultural land. While impact per unit of 
output energy would have been a meaningful indicator to 
compare energy crops, this is not applicable when energy 
and food crops are being compared. The average results (see 
Table 2) are highly variable within a watershed due to the 
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Table 2 Life cycle inventory and life cycle impact assessment results (per square meter per year) of land occupation in ten watersheds of Spain 



Internal watersheds 
—Catalonia 

Ebro 

Duero 

Jucar 

Tajo 

Guadiana 

Segura 

Guadalquivir 

Mediterranean— 
Andalusia 

Atlantic— 
Andalusia 

LCI—soil erosion (10 3 g) 










B-W-R 

1.58 

0.43 

0.27 

1.13 

0.56 

0.31 

0.39 

0.68 

0.69 

0.59 

B-W-P 

1.63 

0.45 

0.29 

1.16 

0.58 

0.31 

0.39 

0.68 

0.68 

0.58 

B-W-F 

2.44 

0.70 

0.46 

1.78 

0.90 

0.50 

0.60 

1.07 

1.01 

0.85 

B-W-OR a 

1.51 

0.41 

0.24 

1.03 

0.54 

0.28 

0.38 

0.63 

0.63 

0.56 

PP a -PP a -PP a 

0.28 

0.09 

0.05 

0.20 

0.10 

0.06 

0.07 

0.11 

0.09 

0.08 

Soil resource 

9 1 

depletion (MJ se m year ) 









B-W-R 

2.5E+04 

6.4E+03 

4.1E+03 

1.9E+04 

9.0E+03 

5.0E+03 

5.8E+03 

9.9E+03 

9.7E+03 

9.9E+03 

B-W-P 

2.6E+04 

6.6E+03 

4.4E+03 

2.0E+04 

9.3E+03 

5.2E+03 

5.8E+03 

9.8E+03 

9.7E+03 

9.9E+03 

B-W-F 

3.9E+04 

1.0E+04 

7.0E+03 

3.0E+04 

1.4E+04 

8.2E+03 

9.2E+03 

1.6E+04 

1.5E+04 

1.4E+04 

B-W-OR a 

2.4E+04 

6.0E+03 

3.6E+03 

1.8E+04 

8.7E+03 

4.6E+03 

5.6E+03 

9.2E+03 

9.2E+03 

9.4E+03 

PP a -PP a -PP a 

4.5E+03 

1.3E+03 

7.5E+02 

3.4E+03 

1.6E+03 

9.0E+02 

9.7E+02 

1.6E+03 

1.4E+03 

1.3E+03 

_o _ 

Ecosystem quality (NPPD m year 

_1 , decimal % of NPP 

o lost m 2 

year -1 ) 






B-W-R 

0.12 

0.03 

0.02 

0.07 

0.04 

0.03 

0.02 

0.05 

0.04 

0.05 

B-W-P 

0.12 

0.04 

0.02 

0.07 

0.04 

0.03 

0.02 

0.05 

0.04 

0.05 

B-W-F 

0.18 

0.05 

0.03 

0.11 

0.06 

0.04 

0.03 

0.07 

0.06 

0.07 

B-W-OR a 

0.11 

0.03 

0.02 

0.07 

0.04 

0.02 

0.02 

0.05 

0.04 

0.05 

PP a -PP a -PP a 

0.03 

0.02 

0.01 

0.02 

0.01 

0.01 

0.01 

0.02 

0.02 

0.02 


B-W-R winter barley-winter wheat-rye, B-W-P winter barley-winter wheat-pea, B-W-F winter barley-winter wheat-unseeded fallow, B-W-OR 
winter barley-winter wheat-oilseed rape, PP-PP-PP poplar-poplar-poplar 

a Crops for energy use 


disparity of soils, climates, and ecosystem biomass productiv¬ 
ities (see ESM Section 6SI). There is a need for statistical 
analyses combining watersheds and geospatial features to 
show to what extent the average results can be extrapolated 
across the watershed. Aggregation at the watershed level was 
not specific enough to reflect a common trend, though it is a 
useful reference area to compare to and/or combine with the 
results of a water use impact assessment. Further highlights of 
the results are discussed in ESM Section 7. 

4.2 Soil erosion impact assessment model 

Analysis of the method followed the general evaluation crite¬ 
ria and specific sub-criteria relevant for land use impacts in the 
ILCD Handbook (JRC 2010). The aim of the analysis was to 
qualitatively address model uncertainty, facilitate comparison 
with other soil erosion impact assessment models, and identify 
the strengths and weaknesses of our method. 

4.2.1 Completeness of scope 

Two relevant impact pathways leading to AoP natural resour¬ 
ces and natural environment were addressed at the endpoint 
level. We expressed ecosystem damages using a growth-based 
value (NPPD), while the majority of endpoint methods express 
ecosystem impacts as PDF (Goedkoop and Spriensma 2001; 


Goedkoop et al. 2009; Jolliet et al. 2003). A significant corre¬ 
lation has been found (Pfister et al. 2009) between vascular 
plant species biodiversity and NPP of the actual vegetation, 
which led the authors to select NPP as a proxy for ecosystem 
quality. The same proxy could be applied to transform NPPD 
to PDF, taking into consideration that we used potential instead 
of actual vegetation. Globally, the NPP 0 -to-NPP act ratio ranges 
from 0.9 to 1.2 on 81 % of the terrestrial area. 

Similarly, for the indicator of soil resource depletion, we 
used emergy units, as in the work by Rugani et al. (2011), 
whereas most endpoint resource demand indicators are 
expressed in energy units (in megajoule), as in the EI99 
(Goedkoop and Spriensma 2001) and IMPACT2002+ 
(Jolliet et al. 2003) methods. In reality, any energy measure 
for reconstructing a renewable resource such as soil lacks 
meaning. To allow for comparisons with other methods and 
impact categories, energy values can be transformed into 
emergies using the resource-specific solar equivalent factor. 
When this is done, the results can be compared to or aggre¬ 
gated with surplus emergy demands of energy carriers (e.g., 
crude oil, 0.091 MJ se g ! ) as well as non-energetic resources 
(e.g., water supply, 0.203 MJ se g -1 ; Zhang et al. 2010). We 
assumed an SEF soil of 23.9 MJ se g -1 soil loss (Odum 1996) 
as this is the only available estimate. This global value does 
not distinguish between different soil types, land uses, and 
world regions, so it should be further refined. 
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One main soil erosion-related impact pathway for human 
health was identified (impact pathway 3 in Fig. 1). This is 
discussed, although not modeled, in ESM Section 8. 

The characterization model and factors are globally ap¬ 
plicable and spatially defined, taking into account ecosys¬ 
tem biomass productivity and soil characteristics at the grid 
level. How to aggregate these grid cell-specific factors on a 
wider and still accurate scale (e.g., ecoregions, land cover, 
watersheds) is a complex unresolved issue due to the huge 
variability of soil types even at the landscape scale. This 
variability of soil types makes it difficult to identify a 
standardized approach to address spatial differentiation. 

4.2.2 Environmental relevance 

The environmental impacts of soil erosion are manifold and 
affect both human life and ecosystems; hence, soil erosion 
control is an important ecosystem service that should be 
considered in the evaluation of the environmental sustain¬ 
ability of land use activities. The method we developed 
focuses on the assessment of land occupation impacts (i.e., 
the use of a land area for a specific human purpose) on soil 
loss for any type of human activity, whenever the soil 
erosion rate of the assessed land use activity (e.g., industrial 
and mining) is available or estimated by the LCA practitioner. 
However, the method discounts impacts due to land transfor¬ 
mation (i.e., change of a land area to make it suitable for a 
specific use). For occupation, it is debatable what period of 
occupation should be considered in the agricultural LCA 
inventory (i.e., duration of the crop, duration of the crop plus 
the fallow period, and crop rotation). For both transformation 
and occupation, the choice of the reference situation to mea¬ 
sure the magnitude of the change and the time needed after 
occupation to recover this reference situation are two funda¬ 
mental issues to be agreed upon to properly assess impacts 
from both land use interventions. So far, the recommended 
reference situation in attributional LCA is the so-called poten¬ 
tial natural vegetation after land occupation (Mila i Canals et 
al. 2007a), which differs from the natural situation because 
nature rarely returns to its original state after being disturbed. 

Of the overall cause-effect chain for land use (see Fig. 1), 
we focused on the soil erosion impact mechanism by quan¬ 
tifying changes in topsoil preservation and changes in the 
NPP 0 due to altered soil function. For a complete evaluation 
of land use impacts, the soil erosion assessment must be 
complemented with indicators that measure effects on cli¬ 
mate change, biodiversity loss, water and nutrient regula¬ 
tion, and unique landscapes. 

4.2.3 Scientific robustness and certainty 

The two indicators reflect the cause-effect chain from inter¬ 
ventions to the final environmental damages. For soil 


resource depletion, the impact can be expressed at the mid¬ 
point level if soil losses of the LCI are weighted by the 
available reference-corrected soil depth without further 
transformation to emergy units. The cause-effect chain for 
ecosystem quality directly links soil erosion to the NPP 0 
depletion indicator, which can be further modeled to PDFs. 

The geographical differentiation of the model has good 
potential for being improved and further developed when 
more detailed global maps of soil properties and biomass 
productivity of ecosystems become available. The current 
large-scale maps do not fully resolve the real diversity of 
soil and ecosystems. 

The indicator for ecosystem quality can be partially veri¬ 
fied against monitoring data by measuring the loss of NPP at 
different degrees of soil erosion. Note that we used potential 
NPP, so the model can be tested in areas with potential or near 
potential vegetation. The indicator for soil resource depletion 
cannot be verified (MJ se ) as emergy is an abstract concept. 

Model uncertainties were qualitatively evaluated using 
the set of criteria listed in the ILCD Handbook (JRC 
2010). Most statistical and decision rule uncertainties still 
have to be estimated. Considerable statistical uncertainties 
occur in the ecosystem quality impact assessment model. 
Our approach of linking soil erosion to NPP loss is based on 
limited and highly uncertain studies about the effect of soil 
erosion on an ecosystem’s biomass production. This adds 
substantial uncertainty to our model. Moreover, we allocat¬ 
ed quantitative ranges of soil loss to qualitative classes of 
soil loss (ESM Table 2SI) and then calculated a linear 
regression with the average of each range. Choosing other 
quantitative ranges or a value different from the average of 
each range would have given different NPPD equations. 
Another source of uncertainty arises from the type of rela¬ 
tionship established between SOC losses and NPP 0 losses. 
We assumed that both variables are linearly related. The 
linear model was used instead of nonlinear models because 
there is insufficient knowledge of the type of relation be¬ 
tween the two variables and the possible interference of 
other variables (e.g., climate change). 

Further statistical uncertainty that affects input data in the 
LCI and the characterization factors is the uncertainty aris¬ 
ing from the resolution of soil data. Also, for the soil 
resource depletion indicator, the use of a site-dependent 
solar energy factor would reduce uncertainty. Uncertainty 
also comes from the type of spatial aggregation of the LCI 
and characterization factor results. These are key future 
research directions to improve the reliability of the model. 

4.2.4 Documentation, transparency, and reproducibility 

The documentation used for the model is published and read¬ 
ily accessible. The maps used to derive the characterization 
factors are available online on valid and reliable sources of 
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data, namely, the FAO’s site (FAO/UNESCO 2007) and in 
published, peer-reviewed, scientific articles (Haberl et al. 
2007). This availability of input data allows third parties to 
further develop and improve the impact factors and the model. 

4.2.5 Applicability 

At the present state of development, two obstacles hinder 
the applicability of the method across the overall life cycle 
of a product. At the inventory stage, unit processes of 
existing LCI datasets currently lack soil loss figures. This 
makes it difficult to include the background system in the 
soil erosion assessment. At the impact assessment stage, we 
regionalized CFs on a gridded spatial scale, while CFs 
should be aggregated to cope with data gaps of location of 
processes. To overcome these problems, processes of the 
LCI datasets should include information on soil erosion, and 
characterization factors, currently disaggregated at the grid 
cell scale, should be aggregated on larger scales. A site- 
generic CF, which is useful for processes of unknown loca¬ 
tion, may be derived by aggregating all the gridded-based 
CFs. The freshwater ecoregion regionalization approach, as 
suggested in Koellner et al. (2012), might be an option of 
aggregation to be studied in future work. 

The CFs are applicable by general LCA practitioners. 
Incorporation in current LCA software will require adaptation 
of the software to accommodate spatially explicit data in GIS 
format. Estimates of the soil erosion data needed for LCI and 
measured with USLE require some knowledge of soil sciences 
and experience in the application of the equation. 

4.2.6 Stakeholder acceptance 

For the ecosystem quality indicator, the results are easily 
interpretable and understandable by non-LCA experts. In 
contrast, a possible barrier hindering acceptance of the soil 
resource depletion indicator is the use of surplus emergy 
units, which may be deemed too complex and meaningless 
for most people. The absence of soil and climate data 
availability can also hamper the applicability and acceptance 
of the two indicators. 

5 Conclusions 

We developed a globally applicable, spatially differentiated 
LCIA method to account for land occupation impacts in 
LCA, focusing on the aspect of soil erosion. The LCI data 
required and data sources and models that can be used to 
obtain the inventory flows were also identified. Spatially 
explicit damage factors on a grid cell-level resolution (ap- 
proximately 10 x 10 km ) for the entire world were provided 
for soil resource depletion and ecosystem quality endpoints. 


The LCA model was successfully applied to agricultural 
plots in Spain to compare soil erosion-related environmental 
impacts that may result from substituting traditional food for 
energy crop rotations. Results from the case study show that 
the lowest erosion rates and environmental damage occur 
when rotating poplar short rotation coppice in Northern and 
Central watersheds in Spain. 

Further research should focus on applying the method 
across the overall life cycle of a product and to identify a 
feasible and relevant spatial scale at which to aggregate 
characterization factors to cope with data gaps on the loca¬ 
tion of processes. Other areas for future research deal with 
improving the reliability of the model by quantifying and, 
insofar as it is possible, reducing statistical uncertainty aris¬ 
ing from the key sources discussed in the article. 
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